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1.0 SUfftARY 

This design note describes the sequential probability ratio test 
(SPRT) for 2-IMU FDI. The SPRT is a statistical technique for detect- 
ing and Isolating soft IflU failures, originally developed for the 
strapdown Inertial reference unit (SIRU) and later adapted at Draper 
Labs to redundant glrballed I?1U's. This note docun'ents the equations 
and describe:, the version of SPRT that will be used for analysis at 
OSC/MDTSCO. Though some of the theory is discussed, a complete 
coverage of SPRT theory is beyond the scope of this report and Is 
relegated to the reference material cited. The flov«hart of a sub- 
routine Incorporating the 2-IMU SPRT Is included, and Is referred 
to in the text for Illustration purposes. Neither test case data 
nor performance evaluation is included, as these will be published 
separately. 

2.0- INTRODUCTION 

Last Septenter the SPRT algorithm v’as baselined at the Level R 
OFT Entry SDR to perform the onboard 2 and 3 IMU FDI testinn diirino 
shuctle entry. In order to both develon and verify the method, a 
SPRT subroutine v/as added to the IMUFDI triple - IMU simulation 
program on the JSC Uni vac 1110 computer. The version of SPRT documented 
In this report has evolved out of the development v/ork to date. Future 
developments in SPRT theory are anticipated and v/ill be documented in 
forthcoming reports. The objective of the studies to date has been to 
optimize performance and sensitivity of the algorithm. With this approach, 
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the narqins between actual perfomance and external performance 
requirements may be easily established. 

3.0 !MPLEMEf!TATIOri DESCRIPTIOri 

This section describes the SPRT alnorithm for nerfomino the 
2-IfiU FDI function with skewed shuttle IWI's. The SPRT subroutine 
of the IHUFDI program (flovfcharted in Appendix 1) is one of several 
versions of the SPRT algoritta that have been explored in the past 
two years. The version described in this report is felt to he both 
the most theoretically straightforward version, and also the version 
most akin to the data tracking test (References 1 and 2), thereby 
enabling an accurate comparison botv-^een the two methods. The Draoer 
Labs version (Reference 3) differs slightly from this implementation 
in some areas. 

The basic idea of this FDI method is to geometrically resolve 
a failure direction to one of four uniouely oriented accelerometers 
or gyros contained in a pair of skewed IflU's (Reference ^). Each 
IMU has one planar (XY) and one sinole axis (7) accelerometer, for 
a total of 4 accelerometer instruments per T’U pair. In addition, 
each IMU has 2 planar gyros (XY, ZR) v/ith the 4th redundant (R) 
axis oriented In the plane of the XY gyro, 12“ from the X axis. 

It is assumed for the 2-IIVi SPRT algorithm that the 7R gyro can 
fail only along the Z axis, since any failure with a component 
along the R axis would be detected by the redundant gyro monitor test 
(Reference 5). Gyro data from the 7 axis can therefore be treated 
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as if it were sensed by a single axis gyro. This assumption brings 
the gyro set into correspondence with the accelerometer set (i.e., 
one XY and one Z axis instrument per IfU), thereby allowing the 
same algorithm to be used for both accelerometer and gyro FPI. 

Ther? are minor differences in error measurement calculation and 
thresholding v/hich will be described subsequently. 

Figure 1 illustrates the functional flow of the 2-IfH' SPRT. 

The method is applicable to any pair of IfHJ's; hovjever, for the 
sake of notational clarity, IMU*s fl and ‘2 are singled out in this 
report. Major functions are discussed individually as follows: 

3.1 Error tleasurenent Calculation 

This function takes output data from the IMU's (tV's,girbal 
angles) and forms a vector equal to the discrepancy (or error) 
bet\ieen two I’tU's. This error vector is then output in the coordinate 
system of each IMU stable platform. 

3.1.1 Gyro 

The gyro error measurement is the "total relative misal innment" 

vector, consisting of the off-diagonal elements of the matrices Q 

and C as shov/n in Figure 2. a. The variables referred to in this 

figure are described as follows: 

B. (n) = is the Euler angle transformation matrix from TRI 
JfJ 

stable platform ?j to the navination base, computed 
using the set of gimbal angles read from IflU at 
the n^^ time step. 
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Figure 1. 2-IMU SPRT Block Flowchart 
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T 

&12 (^) * flimbal angle derived 

transformation from IMU #1 to IWJ #2 stable 

platforms at the n^^ time step. 

Q = ^12^") is the skew symnctric matrix v»hose 

off-diaqonal elements form the total relative 

misalignment vector (fl ) in ?fm *1 coordinates. 
^ u,v,w 

C = is the skew symetric matrix whose 

off-diagonal elements form the total relative 
misalignment vector (O ) in IMM ?2 coordinates. 

X »y jZ 

In the real time computer the six small misalignment angles, 9u,v,w 

and 0 , can be computed directly as dot products of rows and 

x»y>z 

columns of the transformations 

3.1.2 Accel erometer 

The accelerometer error measurement is the "incremental AAV" 
vector formed by differencing IMU sensed AV's in the coordinate systems 
of each IHU cluster, as shown in Figure 2.b. The current oimbal 
angle derived transformation, is used for all transformations 

to minimize the effects of oyro drift on accelerometer FPT. 

3.2 V/hiteninq Filter 

A first order recursive filter is used to transform the correlated 
IMU error measurements in to a sequence of independent samples. The 
theory for using a filter in this v/ay is enhodied in the widely knov/n 
result from estimation theory that the residuals of an optimal Kalman 













-7- 


f 11 ter constitute a white Gaussian noise sequence (References 6 and 
7). In this SPRT implementation, each axis is filtered separately 
by a single state estimator. The outputs of the filters are the 
residuals formed by subtracting the filter estimate from the actual 
IMU error measurement. The v/hitening transformation performed by 
this filter is required because independancy of samples is assumed 
in formulating the LLR update equation (Section 3.^). 


3.2.1 Filter Equations 


First pass initialization: 

^ .g(-DELTAT/T) 

State propagation 

^n ' ^ ^n-1 
Residual computation 


r = Y - S 
n n n 


State update 


= S„ + K r„ 
n n n 


.th 


where S„ = state estimate at n time step 
n 

r^ = residual at n^^ time step (filter output) 

= IMU error measurement at n^” time step (filter input) 

K = filter gain (pre-mission c istant) 

T = Autocorrelation time constant (pre-mission constant) 

The above filter is used for both gyro and accelerometer data; the 
only difference being in different state update gains (K), v/hich are 
listed in section 4.2. 
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3.2.2 Pre-mission Filter Tuning 

The filter gains mentioned above v/cre computed pre-mission 
via a tuning procedure derived from Mehra, Reference 6, which 
basically states that the optimality of the filter can be measured 
by how uncorrelated the residuals are. In a routine v/ritten for 
the Univac 1110 Demand terminal, the filter gains v;ere varied 
until the average autocovariance of the filter residuals over 
30 Monte Carlo cycles v/as zeroed out. The procedure used for 
calculating the average autocovariance follovfs: 


Step 1. Calculate a residual mean for each of f axes in each 
Monte Carlo Cycle between the times 3^5 and 12<’5 (this 
time period contains no vehicle attitude transients): 


1 


1295 ^ 

^’axis, cycle “ TM ^ 

bv 5 


Step 2 . Calculate the first autocovariance coefficient, P.(l), 
for each of 6 axes in each Monte Carlo cycle: 

1295 

JL T. (r.-b) (r* c-b) 

198 t=310 ^ ^ ^ 

by 5 


^axis, cycle 1 « 

F99 5 : (r.-b)^ 

t=3H6 ^ 

by 5 

Step 3. Calculate the average autocovariance coefficient. 


ave 


= i_ E 

TCO axis=l 


30 

E 

cycle=l 


^axis, cycle 




-9- 


Repeat steps 1 through 3 on residuals acquired by iterating on the 
filter gain K until a of approximately zero is obtained. 

V/hile computing the gyro gain of section 4.2 = *54) 

the following autocovariances v/ere obtained: 

^ave " 

^min (X axis. IMU cycle 420) 

Pn,,v •= 64.2?: (X axis, VV 42, cvcle «27) 

HiaA 

While computing the accelerometer gain of section 4.2 ~ 

the followina autocovariances v;ere obtained: 

Pave = 

Pmin " -62.0^ U axis, IMU 42, cycle 423) 

Pr«v = 81.9% (Z axis, IMU 42, cycle 430) 

iIIq A 


3 . 3 ILR Parameter Calculation 

The LLR update described in section 3.4 requires tv/o parareters 
to compute the LLR - the classification threshold and the standard 
deviation for the filter residuals. 


3.3.1 Residual Threshold 

The SPRT algorithm requires thresholds separating residuals 
indicative of failure operation from residuals encountered durinn 
nominal operation. As in all threshold type tests, it is desirable 
to have the threshold tiohtly fitted for oreater sensitivity, yet 
still high enouoh so that no false alarms are produced under nominal 
operation. The failure residual thresholds used in this version 
of SPRT were generated by multiplying a base failure threshold 
(Figure 3) by a dynamic scaling factor which accounts for the filter 




ffop l«Ujr li«tcrf«ca 


Figuire 3.b) Incremental AV Differences 


ToucMomi 
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effects of bias attenuation and attitude transients. The base 
thresholds, specified in section 4.4, are the sane closely fitted 
thresholds used in the data trackinn test (Reference 1). The 
scaling factor r » derived fron Reference 8, is calculated by: 


* rn ' rUB.r.) 


where: 


l-<^ 

^ss ‘ i-i(i-ky 


letting n=t/At, 


^n ■ !-<.(! -K) 

K(J) 

%,n " 1-(!>(1-K) 


(*(l-k))"‘J‘ 


in v/hich 

K = filter gain 

PjlB ' P®*^centage attributed to vehicle 
attitude transients (Section ^.5) 
j is set equal to n v/henever the create';! 
ginbal angle change exceeds the limit 

^Gim 

The final threshold is then 


BR(t) = r(t) T(t) 


which is plotted in Figure 3. 


3.3.2 Residual Standard Deviation 

In addition to the residual threshold, the SPRT requires 
a residual standard deviation characterizing the "spi'ead" of 
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nominal data. Standard deviations v/ere conputed at several tine 
slots by the formula 


oit) =/, 


6 30 

18 0 ^ , r. 

axis=l cycle=l 




1 

8 


t-. ; 

E 

t'=t 
by 5 


1/2 


Residual standard deviations \ierc found to be fairly constant at 
the values given in section ^t,5. 


3.4 Loo Likelihood Ratio Update ? 

The likelihood rat’o is the heart of the SPRT alpori'hn, 
being proportional to the probability that a failure has ru red. 

The loo of the likelihood ratio is used to avoid computation of 
the EXP function at each time point. Lon likelihood ratios (Llf's) 
are computed sequentially, modciing the filter residuals as i dep.-ndent 
Gaussian samples (not necessarily zero mean', using the residual 
thresholds and sigmas computed in section 3.3. In this implementation 
the LLR for each axis is split itito tv/o suh-l.LR's, one for positive 
and one for negative failures, since either polarity is assumed 
equally likely. The HR equations follov/: 


LLP. Upda te 


. 2DR 

^ - ^-1 ^ r ? 


(t) 


(r. 


CR(t)) 
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LIR Reset 

If A* < 0 then A* = 0 
n n 

If A~ < 0 then A~ = 0 
n n 

Take the biggest: 

A„ = max (A*, A~) 
n n n 

3.5 LLR Threshold Calculation 

The LLR failure decision threshold is calculated according to 
Reference 9, using the residual thresholds and sigmas from section 
3,3, by the following equation: 

Bf.'(t) = In (2BR^(t)/ALPHA 

3.6 Single Axis/Pual Axis FDI Lonic 

This FOI function performs a threshold test on the f error 
neasurenent channels for each instrument type, both accelerometer 
• and gyro. Each instrument (200F, IDOF) is tested independently 
for evidence of failure. For example, if either the X or the Y 
gyro neasurenent from I!*U *1 is out of tolerance, the failure has 
been detected in the XY^ gyro. Figure 4 shovfs the logic by which 
individual instrum.''nt failure detection test are combined to iso- 
late the failure to a particular instrument. An IHU failure has 
teen detecte ' v/hen at least one of the four instruments of either 
type* or accelerometer, is out of tolerance. In the IfMiFPI 
implementation, at least tv>o simultaneous threshold crossinr? out 
of the four instruments is required before an I*’M failure detection 
is registered. A single threshold crossing v/ould be sufficient. 





Figure 4. Single Axis/Dual Axis FDI Logic 
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but is also potentially mere susceptible to false alarms. A failure 
can be isolated when exactly three of four instriments simultaneously 
show out-of -tolerance mcasuranents. The v»here the failure has 
been detected in both instruments is the unfailed IM!I, since tvfo 
detections vxjuld be indicative of a simultaneous failure - an assumed 
impossibility. Consequently, the If*U in which the failure has been 
detected on only one instrument neasurenent is the failed IMU. The 
case of four simultaneous threshold crossings (not covered in the 
flowchart of Figure 4) is an abnormal condition and should be flagged 
to the crew. 

4.0 RESULTS 

The following paragraphs present a group of constants, pertinent 
to the 2-P1U SPRT algorithm, that vrere used in the r*!iFDI prooram 
to generate the results contained in Reference 1^. These constants 
apply to the reference mission 30 entry trajectory. 

4.1 Time Step 

DELTAT = 5 sec. 

4.2 Filter Constants 
Autocorrelation time constant: 

■^GYRO ^ '^ACCL "" 

Filter gains: 

‘^YRO 
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4.3 St art - Stop Tines 

W'® (''"‘'•o®® «•) 

Tgnd ' (Touchdovm) 

4.4 Base Failure Thresholds 
Gyro: 

3 

TG(t) = X) C. t^ (radians) 
i=o ' 

Cq * 3.7000 E-4 
* 3.4967 E -6 
Cg =-1.1786 E-10 
C3 =-5.8263 F-13 
Accelerometer: 

for 0<t^n45, 

3 

TA(t) = X V (iTi/sec) 
i=o 

Cjj = 2.8300 E-5 
= 4.3503 E-8 
C2 =-5.3665 E-11 
C3 = 9.7743 E-14 
for t>11^5, 

TA(t) = TA(t=1145) 

4 . 5 HR Constants 
Gyro: 

ctq = 2.4 E-4 radians 
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At;cclcrometer: 

PriB * 

=1.2 E-5 Km/sec 

4.6 Mean Tine I’etv/cen False Alams 


T=5000 sec 

lean False Alam Rate: 
ALPHA=DF.LTAT/T=in‘^ 
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APPEUDIX 1 


Flowchart of the SPRT Subroutine 
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SUPROUTiriES CALLFO PY SPRT 


1) o^rrxfi M3) 

Ftetrix transposed tirts a matrix. 

"3 • ”1 ''2 

2) OMxni (M^. M3) 

Matrix tines matrix transposed 

Hj . M, n\ 

3) DMTXV («,, Vg, V 3 ) 

Matrix transposed tines a vector 


V. = k: V, 

o I c 


4) OMXV (M^. V2, Y3) 

Matrix tines a vector 

V = r* V 

'3 1 ^^2 

5) OVSUP (V^ , Y^. Y3) 

Vector subtraction 


The above subroutines assune that all vectors and matrices are 
dimensioned (3) and (3,3), respectively. 

6) Tlf'.ESl and TIHES2 are used by the driver pronran to provide a 
surnary of detection and isolation tines at the end of a rult 
Monte Carlo cycle run. 

7) DVSCLH (SCALAR^, V^, V3) 

Scalar tines a vector 


V3 = SCALAR^ 
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8) DVAOO (V^ , V^, V3) 
Vector Addition 



+ V 


2 



